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Abstract 

In this paper we discuss the issue of geometric stiffening as it arises in the context of 
multibody dynamics. This topic has been treated in a number of previous publications 
in this journal and appears to be a debated subject. The controversy revolves primarily 
around the “correct” methodology for incorporating the stiffening effect into dynamics 
formulations. The main goal of this work is to present the different approaches that have 
been developed for this problem through an in-depth review of several publications dealing 
with this subject. This is done with the goal of contributing to a precise understanding of 
the existing methodologies for modelling the stiffening effects in multibody systems. Thus, 
in presenting the material we attempt to illuminate the key characteristics of the various 
methods as well as as show how they relate to each other. In addition, we offer a number 
of novel insights and clarifying interpretations of these schemes. The paper is completed 
with a general classification and comparison of the different approaches. 


1 Introduction 

The issue of geometric stiffening, also referred to as dynamic stiffening, centrifugal stiffening 
and foreshortening has been a topic of many recent publications dealing with the dynamics of 
flexible bodies for applications to multibody systems. Kane et al. 1 first observed that the 
majority of existing multibody dynamics formulations and accordingly the dynamics simulation 
packages do not incorporate the geometric stiffening effect. They have attributed this flaw to 
the “convenional approach” for describing the deformation of elastic bodies, which yields a set 
of dynamics equations that inherently lack the geometric stiffening terms. Kane et al. proposed 
an alternative approach, correcting this flaw, and applied it to develop a set of equations 
for the deformation of a beam attached to a moving base. 

Eke and Laskin 2 took up the issue raised in Ref. 1 and investigated regimes of validity of 
existing formulations with the simulation package DISCOS on a spin-up beam example. They 
qualified the error in conventional approach as a “premature linearization” of the displacement 
field. This was later supported by Padilla and von Flotow 3 and Banerjee and Dickens 4 . 

Shortly after Kane et al.' s publication, two commentaries appeared on the material 
presented in Ref. 1. In particular, London 5 pointed out that geometric stiffening has been 
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previously considered by many researchers in a number of applications. Furthermore, he 
observed that several approaches have been employed to include this effect in the dynamics 
equations and compiled a table characterizing the various methods. In a technical note , 6 
Hanagud and Sarkar state that, contrary to the claim made in Ref. 1, the conventional method 
for modelling the kinematics of elastic deformation can be used. 

The existing controversy over the nature of geometric stiffening, the debate on “ the 
correct” approach to model it and the seeming incongruity of the methods used to include the 
effect in the motion equations — all of these have motivated us to review several of the works 
on this subject. In doing so, we have attempted to understand precisely how geometric stiff- 
ening is incorporated into the dynamics equations in different approaches, what assumptions 
and approximations are made in the derivation, what motivated these and whether they are 
justified. This paper contains the main results of our review. 

Our starting point will be the landmark paper by Kane et al. 1 and the subsequent 
commentaries . 5,6 Following that, we give a thorough treatment of the works by Likins et al., 7 
Vigneron 8 and Kaza and Kvaternik 9 and a summary of the relevant material from the pub- 
lications by Lips and Modi 10 and Hughes and Fung . 11 Section 4 contains the main results 
from Laskin et al., 17 Meirovitch 13,14 and Banerjee et al ,. 4,15 In reviewing the works of these 
researchers, we do not simply repeat their derivations, nor do we include the dynamics equa- 
tions developed in these publications. Instead, we concentrate on the fundamental assumptions 
made in formulating the basic elements necessary for deriving these equations, where the “for- 
mulation” ends when the development becomes a purely mechanical process. For instance, in 
the cases where dynamics equations are derived via Hamilton’s principle, we limit ourselves 
to stating kinetic and potential energy functions, and do not go through the procedure of 
applying the variational principle. This allows us to compare the various approaches based on 
the fundamental physical assumptions. 

In addition to presenting the key features of different procedures, making comparisons 
and establishing relationships between them, we also provide clarifications and give some new 
insights. We conclude the paper with a discussion in which we disclose some of the exist- 
ing misconceptions, classify the approaches and comment on their suitability for multibody 
dynamics simulation. 


2 Kane et al. and Commentaries 

2.1 Main Results of Kane, Ryan and Banerjee 

In Ref. 1 , Kane, Ryan and Banerjee develop the dynamics equations of a general flexible beam 
built into a rigid base. The base body can undergo arbitrary, but prescribed translational 
and rotational motion. The generality of the beam refers to the fact that its geometric and 
material properties are not assumed to be constant, but can vary along the length of the beam. 
In addition, Kane et al. do not make the common assumption that the elastic and centroidal 
axes coincide. As a result, their motion equations contain terms dependent on the components 
of the eccentricity vector, e-i and e 3 . 

The formulation of equations in Ref. 1 differs from many existing procedures in several re- 
spects. First, it incorporates the effect of the transverse displacement on the axial displacement 
in the kinematic description of the deformation. This is achieved indirectly by expressing the 
distance along the deformed elastic axis as a nonlinear function of the transverse displacements 
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with: 


H£)’ *(fc)l 


1/2 


da 


( 1 ) 


The above equation is the same as Eq. (19) of Ref. 1, where the variable s denotes “the stretch 
in the beam along the elastic axis.” 

The second major difference relates to the choice of elastic deformations that are em- 
ployed to describe kinematics of the deformed beam. The standard procedure is to use three 
orthogonal elastic displacements ui, u 2 and u 3 to represent the displacement field in a de- 
formable body. Kane et al employ the stretch s with the transverse translations u 2 and u 3 
as a set of generalized coordinates. Thus, when discretizing the continuous displacement field, 
they discretize the stretch variable instead of the axial displacement uj. This is expressed by 
Eqs. (25) and (26) of Ref. 1, which we rewrite here for convenience as: 


s(x ,t) = 

j = l 

Ui(x,t) - i = 2,3 

j=l 

( 2 ) 


Accordingly, the dynamics equations based on the above premise reprsent a model for the time- 
evolution of {s, u 2 , u 3 }, or rather, the corresponding discrete elastic coordinates. The “conven- 
tional approach” involves discretizing the orthogonal set of elastic displacements {ui,u 2 ,u 3 } 
with: 

u t = J24> t] q } , * = 1,2,3 (3) 

i=l 


or in matrix form: 


Ml 

u 2 

«3 


= u = 


( 4 ) 


Kane et al. argue that in the standard procedure (3), the three elastic deformations cannot 
account for the fact that every transverse displacement gives rise to an axial displacement, 
because the form (3) inherently precludes such an interdependence. 

The general methodology employed in Ref. 1 to derive an explicit (literal) set of motion 
equations for the elastic coordinates is that presented in Kane and Levinson. ^ The procedure 
requires one to construct the generalized inertia and generalized active forces. The former are 
developed in Ref. 1 according to the algorithm outlined by Kane and Levinson. The generalized 
active forces, which for the particular system considered result from internal forces, are derived 
from the strain energy function. The expression for this function U is given in terms of the 
components of the force and torque vectors which act on a cross-section of the beam (see Eq. 
(51) in Ref. 1). Thus, it takes the form of a sum of six integrals, three for each of the force and 
torque, where the integrand of each integral term is a quadratic function of the appropriate 
load. In order to determine the generalized internal force by using Eq. (50) of Ref. 1, which 
in fact is a statement of Castigliano’s theorem, one needs to formulate the strain energy as a 
function of the generalized coordinates. To this end, Kane et al. express each of the six loads 
as a linear function of elastic deformations or their spatial derivatives. 
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For reasons that will become apparent in §2.4, we draw attention to one particular term 
in the strain energy function. This term represents the contribution of the axial load and will 
be denoted here by Up. It corresponds to the first term in Eq. (51) of Ref. 1, which we rewrite, 
omitting the subscripts as: 

f L P 2 

U '=!o2EA ix (5) 

In linear analysis, one approximates the axial load as a linear function of the axial displacement 
gradient with P = EA^£. The expression for P given in Ref. 1 by the left Eq. (58) takes the 
same form, but with the axial translation u\ replaced with the stretch variable s. With that, 
the axial load contribution to the strain energy function becomes: 



( 6 ) 


The above are what we view as the key features of the dynamics formulation for a flexible 
beam attached to a moving base that has been put forward by Kane, Ryan and Banerjee. 1 In 
the following two subsections, we summarize the main points of a technical comment and an 
engineering note, both of which are related to Ref. 1. These appeared in two issues of the 1988 
volume of the Journal of Guidance and Control, shortly after the publication of Kane et ai 
Section 2 is concluded with a discussion of the two commentaries. 


2.2 London’s Comments 

In summarizing London’s comments, we have grouped them into two categories. The first one 
includes comments which deal with the “qualitative” aspects of Kane, Ryan and Banerjee ’s 
work, such as literature review. In the second category, we include comments related to the 
quantitave or technical aspects of the analytical development presented in Ref. L 

Category I. London observes that Kane et ai “create the impression that a new theory has 
been discovered,” refering to the theory to model foreshortening. London bases this 
statement on the fact that Kane et ai do not give any references as to the origin of their 
results for the nonlinear description of the kinematics of the deformed beam (primarily, 
Eq. (19) in Ref. 1). 

Regarding the literature review on the subject of modelling flexible beams attached to 
a moving base, London’s criticisms are twofold. First, he suggests that the references 
in Kane et a/.’s manuscript are incomplete. To be specific, they neglect to mention the 
work by Lips and Modi 19 on the dynamics of beams undergoing a three- axis spin, as 
well as the work by Hughes and Fung, 11 which treats the problem of stability of spinning 
satellites with flexible appendages. Second, London implies that many of the references 
that are included have not been given a proper and/or appropriate credit. In this regard, 
he particularly notes the work of Kaza and Kvaternik 9 which is classified by Kane et 
ai into a group of papers “in connection with aircraft dynamics” addressing “questions 
concerning tapered, twisted and rotating beams.” Similar treatment is given to the works 
of Likins et ai 7 and Vigneron, 8 which are grouped under those in the field of “spacecraft 
dynamics” with a “particular interest” in “the effect of vehicle elasticity on attitude 
motions.” 

Category II. With regards to the technical merits of Kane et a/.’s formulation, London draws 
attention to the following four points: 
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1. London questions why the authors of Ref. 1 “choose to represent the three elastic 
degrees of freedom by five variables {uj, u 2 ,U3, s,£}.” 

2. It is pointed out that the use of x + U\ as an upper limit of the integral in Eq. (1) 
(Eq. (19) in Ref. 1) creates problems in evaluating the modal integrals /3,j and 7,j 
(Eq. (41) in Ref. 1). 

3. London questions why the foreshortening is not considered during evaluation of the 
strain energy. 

4. It is stated that the final equations derived by Kane “are linear in terms of vibration 
coordinates (but still retain higher-order spin effects) as is already the case with most 
other works.” 


2.3 Hanagud and Sarkar’s Note 


Contrary to what is claimed in Ref. 1, Hanagud and Sarkar believe that the axial and transverse 
motions can be treated independently with the standard discretization procedure and the 
stiffening effect can be accounted for if “the nonlinear effects are properly included in the 
formulation.” 6 

Hanagud and Sarkar present a formulation where, as in the conventional approach, they 
discretize the axial displacement u u not the stretch variable s. Hanagud and Sarkar derive 
the differential equations for the corresponding discrete elastic coordinates by using the same 
general methodology as employed in Ref. 1 , which as they point out is sometimes referred to as 
Kane’s method. Similar to the development of Kane et al . , they also determine the generalized 
active forceg which are the elastic (internal) forces, from the strain energy function. However, 
Hanagud and Sarkar formulate the strain energy as a quartic function of the spatial derivatives 
of ui, u 2 and U3. This is accomplished by employing the nonlinear strain-displacement rela- 
tions, through which the aforementioned nonlinear effects are introduced into the formulation. 
Lastly, Hanagud and Sarkar do not linearize the final equations of motion, but retain terms of 
second and third order. 

An important contribution of Hanagud and Sarkar’s work is an observation that the 
expression for the stretch variable presented by Kane et al. (Eq. (19) in Ref. 1) is inconsistent 
with the rest of their development. The inconsistency results from the fact that this relation 
for s is applicable if one expresses the transverse displacements u 2 and U3 as a function of the 
deformed coordinate X which corresponds to the axial projection of a generic point in the de- 
formed configuration of the beam. (We have chosen to follow the traditional notation employed 
in the theory of elasticity, where one distinguishes the deformed and undeformed coordinates 
by different- case letters.) Since in their formulation, Kane et al. express the translations u 2 
and U3 in terms of the undeformed coordinate x (see Eq. (2)), the consistent expression for the 
stretch variable is: 



dui(x,t) 

dx 



duj (s,Q 

dx 



du 3 (x,t) 

dx 



21 1/2 


do (7) 


The above equation, although looks different, is equivalent to Eq. (2) in Ref. 6. We also 
observe that it embodies the nonlinear formulation of the strain. In (7), we have used the 

notation to emphasize that elastic translations must be expressed as a function of the 

undeformed coordinate x, while a in this case is a dummy integration variable. We also point 
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out "a change in the upper limit of the integral from x + U\ in Eq. (1) to x in Eq. (7) which, of 
course, is a consequence of the transformation from deformed to undeformed axial coordinate. 


2.4 Discussion of the Commentaries 

In the following, we offer our opinions on the various points made in the commentaries. We 
hope that these will serve to elucidate the more subtle features of the formulation proposed by 
Kane et ai A couple of the issues raised in the commentaries will be further addressed in the 
final section of the paper. 

We agree with London in saying that when one reads the manuscript, 1 one gets an impres- 
sion that its authors propose a new theory to model the deformation of a beam. Fundamentally, 
the theory is not new and a number of formulations which incorporate the foreshortening effect 
and the resultant stiffening of the beam during Its rotation, have been previously published. 
(The material presented in §3 and §4 will support this fact.) Indeed, Kane et ai refer to 
some of these works, but only in a superficial and in some cases misleading manner. This 
notwithstanding, we feel that the formulation developed in Ref. 1 does have a couple of novel 
features. These are: (i) generality of the system modelled and (ii) use of the stretch variable as 
a generalized coordinate in deriving the dynamics equations. Although we do not share Kane’s 
conviction that the dynamics equations must be formulated in terms of the stretch variable in 
order to predict stiffening of the beam, this particular feature of their procedure provides, at 
the least, an interesting alternative to the conventional approach. 

Continuing with London’s comments in the second category, we offer the following ob- 
servations: 

1. In our interpretation of the formulation, 1 the stretch s is introduced to replace the axial 
displacement u\. This is made abundantly clear in section IV of Ref. 1 where s takes 
place of what usually appears as «j. The variable £ is employed as a short-hand for the 
combination x + u\. In this light, we do not agree with London’s statement that Kane 
et ai propose to use five variables to represent a three-degree-of-freedom displacement 
field. 

2. According to the development presented in Ref. 1, the modal integrals /? tJ , 7 are by 
definition time-dependent through u x = u^x^in their upper limit of integration. This 
implies that they must be evaluated at each time step in the numerical integration of 
the motion equations. (Although that would certainly add to the computational cost 
of the simulation, it should not pose a problem otherwise.) The situation changes, 
however, if one corrects formulation in Ref. 1 in accordance with observations made 
by Hanagud and Sarkar. That can be achieved by replacing the inconsistent expression 
for the stretch, given by Eq. (1), with a consistent form of Eq. (7). With this correction, 
the aforementioned modal integrals become independent of time and, therefore, can be 
evaluated prior to the numerical integration of the motion equations. 

3. This particular comment by London has motivated us to consider closely the expression 
for the strain energy function employed in Ref. 1. In the process, we have singled out 
the “axial” strain energy, previously denoted by f/p, as the only contribution which may 
possibly comprise the foreshortening effect. The following brief development shows that 
indeed it does. 

The function Up, as can be seen from (6), is quadratic in the spatial derivative of the 
stretch, Thus, in order to illuminate the nature of this strain energy, we need to obtain 
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an explicit expression for the stretch gradient. As noted previously, Eq. (7) provides the 
proper form for s that should be employed in Ref. 1. Abbreviating the notation, we 
rewrite (7) for s with: 




— x 


Differentiating the above with respect to x we get: 


ds _ 

[/, duA 2 , 

fdu 2 \ 2 . 

/ 0 « 3\ 2 

dx 

( ,+ «») + 

05 1 
+ 

[dx ) 


f|i/2 


- 1 


and upon expansion of the first term, the required gradient takes the form: 

1 1/2 

- 1 


ds 

, n dui 

( du\ \ 2 , 

(du 2 \ 2 

( du 3 \ 2 

dx 

1 + 2 -^- + 
dx 

Ur) + 

Ur) + 

[dx ) 


( 8 ) 


(9) 


( 10 ) 


Let us now introduce the axial strain £ 0 ,ir , where the 0 subscript signifies that it refers 
to the elastic axis. (Note that tij , u 2 and u 3 are defined in Ref. 1 as translations of points 
along the elastic axis only.) The strain £o,xx can be expressed in terms of the elastic 
displacements with a well-established strain- displacement relation. It has the following 
exact and nonlinear form: 


£0,xr 



(ii) 


With the above, the stretch gradient of Eq. (10) can be succinctly written as: 

|^ = [l + 2eo,x*] I/2 -l (12) 

Before we continue, it is worthwhile to point out that in all formulations dealing with the 
subject of geometric stiffening in the context of multibody dynamics, it is always assumed, 
although not always stated, that the strains are small, and specifically, £o,xx << 1* 
Therefore, we can make use of the Binomial Theorem to simplify Eq. (12). Retaining 
the first two terms in the binomial expansion we get: 


ds 

dx 


« [1 + £o,xx] - 1 


£o,rr 


(13) 


Finally, substituting for || from (13), the axial contribution to the strain energy function 
employed by Kane et al. takes the form: 



(14) 


with the axial strain given by the nonlinear Eq. (11). 

At this point, it is appropriate to comment on the form of the strain energy employed 
by Hanagud and Sarkar in their formulation (Eq. (8) in Ref. 6). Their expression can be 
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derived by using Eq. (14) in conjunction with the following simplified expression for the 


axial strain: 



(15) 


The above form is obtained from the exact relation (11) by omitting the (7^) term 

an approximation often made in the context of moderate-deformation theories.* 7 The 
third- and fourth-order terms in the strain energy of Hanagud and Sarkar result from the 

nonlinear terms and (^?) in the strain-displacement relation (15). It is these 

terms that lead to the geometric stiffening and foreshortening of the beam. 

The development presented here demonstrates that, contrary to London’s comment, the 
geometric nonlinearity is included in the strain energy expression given in Ref. 1. This is 
not apparent because Kane et al. employ the stretch s instead of the axial translation u\ 
as a generalized coordinate in their formulation. As a consequence, they do not need to 
expand |* when evaluating U so that their strain energy function remains quadratic in 
the discrete generalized coordinates and does not explicitly contain higher-order terms. 

4. We have found this comment by London somewhat bewildering. As stated in it, and as 
clearly stated in Ref. 1, the equations of motion developed by Kane et al. are linear in the 
elastic coordinates and incorporate the geometric stiffening. In fact, one of the crucial 
points made in Ref. 1 is that the linear motion equations in other works do not contain 
all of the linear terms, in particular the stiffening term. As implied by Kane et al., 1 
and later stated by Eke and Laskin 2 and Padilla and von Flotow, 3 this occurs because 
of premature linearization implicit in the “conventional approach.” Contrary to this, we 
support Hanagud and Sarkar’s view that one can obtain the stiffening effect with the 
conventional approach. This is achieved by employing the nonlinear strain-displacement 
relations in constructing the strain energy function. The stiffening term obtained with 
this approach is a nonlinear function of elastic coordinates. 

The reason why the final equations of Kane et al. are linear in elastic coordinates, but yet 
include the stiffening terms lies in their choice of the stretch s as a generalized coordinate 
and the fact that it includes the nonlinear contribution from transverse displacements. 
As we had shown in the previous comment, it is through the use of stretch instead of axial 
displacement, that Kane et al. incorporate foreshortening in their formulation, without 
introducing nonlinearity explicitly into the motion equations. 


3 References Contended by London 

London’s comments on the literature review and Kane et al .' s reply to them 18 , motivated us to 
investigate the material presented in several of the references in question. The main objective 
of this section is to summarize these findings. 

3.1 Likins et al., Vigneron, and Kaza and Kvaternik 

Our choice to discuss the contributions of Likins et al., 7 Vigneron, 8 and Kaza and Kvaternik 9 in 
the same section is based on several reasons. First, these articles appeared within a time-span 
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of four years and therefore belong to the same “era.” Second, all three publications contain 
formulations of the dynamics equations for a flexible beam spinning in a plane at a constant 
speed. Finally, there is a logical relationship between these works, since that by Vigneron is 
a comment on Likins et al., while Kaza and Kvaternik extend Vigneron’s approach to obtain 
the nonlinear equations of motion of the aforementioned system. 

Our presentation is not a plain copy of Refs. 7, 8 and 9, as it is structured to make 
apparent the key features of the approaches taken in the three works, establish a relationship 
between them, as well as identify the particular contributions of each one. Towards this 
end, we present the results of these works in a common notation, which will also be employed 
throughout the rest of the paper. This notation is similar to that used in Ref. 1 with one major 
difference. We choose to denote the three orthogonal elastic displacements with symbols u , v 
and w. Thus, u now represents the axial elastic displacement measured along the x-axis of the 
reference frame, while v and w are the transverse elastic displacements. Furthermore, these 
symbols are not restricted to the elastic axis, but represent elastic displacements of any point 
in the beam. This convention follows that used by Kaza and Kvaternik. 9 

As already mentioned, the system considered in all three publications 1 is a uniform elastic 
beam, with a symmetric cross-section, spinning at a constant angular speed fl about the z-axis 
of the reference frame. This system represents a special case of that treated by Kane et al., 
defined with: u> 3 = ft, = w 2 = = v 2 = v 3 = 0, e 2 = e 3 = 0 as well as constant geometric 

and material properties. The common features of the formulations 7,8,9 are listed below. 

(i) The dynamics equations are constructed via Hamilton’s principle. 

(ii) The position of a generic point located at [x,y, z] T in the undeformed beam is given 
by [x + u,y + v,z + u?] T , where v and w are functions of x and time only, that is 
it(x, y, z, t ) = u(x, t) and similarly for v and w. As well, the transverse displacements are 
assumed constant in a cross-section, thus precluding torsional deformation. 

(iii) The kinetic energy is calculated with: 

T = \ ///[*» + * + * + ft 2 (x + «) 2 + ft 2 (y + v) 2 

+ 2ft(x -(- u)v - 2ft(y + u)iij dx dydz (16) 

where, like in Ref. 1, the symbol p denotes the mass per unit length of the beam and A 
is the cross-sectional area. The above equation is identical to equation (6) of Likins et 
al. and applies to any elastic body spinning as specified before. It can be expanded and 
simplified for a beam with a symmetric cross-section to yield: 

T =\p J (u 2 + v 2 + w 2 ) dx+\fl 2 p J (x 2 + u 2 + v 2 + 2xu) dx 

J (xv + uv — vii) dx+^pfl 2 LI z 

where we have used the standard definition I z = ff y 2 dydz. 

(iv) The potential energy is calculated with 

£/ = f ///(£** )2 d* dydz (18) 

1 Likins et al. also consider “Axial Beams.” 
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where £ xx denotes the axial strain at any point in the beam, and is given by a nonlinear 
expression: 

2 21 
- + 


du 


Sxx ~ dx + * 


(r) + Gr) +(^ 

\dxj \dxj \ox 




(19) 


(v) The dynamics equations are formulated for the continuous displacement variables, and 
accordingly take the form of partial differential equations. 


There are two major differences between formulations presented in Refs. 7, 8 and 9. The 
first one relates to the form of the assumed axial displacement field. Likins et al. “expand” 


the axial displacement u with: 


dv 


U = Uo - y d~x 



( 20 ) 


Vigneron, followed by Kaza and Kvaternik adopt a different form. Their axial displacement is 
given by Eqs. (2) and (la) in Refs. 8 and 9 respectively, which we write as: 


dv 



( 21 ) 


In the above, uj is the “displacement associated with the foreshortening effect”. 8 In both 
references, this component of the axial displacement is specified as an explicit function of the 
transverse displacements: 


Uf = uj(x,t) 



( 22 ) 


We note that expression (22) is a second-order approximation for the foreshortening of the 
beam. Clearly, uj corresponds to the axial displacement which results from the transverse 
deformation. Comparing Eqs. (20) and (21), we observe that «o must equal u a - uj and hence 
« 0 and u a represent different physical quantities. We have employed the subscript s in u 3 to 
signify that it corresponds to axial displacement resulting strictly from extension or stretch of 
the elastic axis. By comparison, uq represents the total axial displacement on the elastic axis 
which may be due to both stretch and bending. 

The second difference between the derivations of Likins et al., Vigneron, and Kaza and 
Kvaternik lies in the approximations made in deriving the final equations of motion for the 
elastic coordinates u 0 or u a , v and w. To be specific, these are approximations made in 
formulating kinetic and strain energies that are subsequently used to construct the Lagrangian 
of the system. 

Likins et al. substitute their expansion for u from Eq. (20) into the kinetic and strain 
energy expressions (17) and (18), the latter combined with (19), thereby reformulating T and 
U in terms of uq, v and w. They simplify the results by making the following assumptions. 


I. The in-plane deformation is ignored based on the argument that it is present only because 
of the Poisson effect. 

(a) This amounts to dropping all terms in kinetic energy T which involve n or its 
derivatives. 
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(b) The strain energy U is derived by omitting in 4h e definition (19) for the 

strain. 

II. Likins et al. neglect all terms that involve u 0 and u 0 in T, thus effectively eliminating 
the axial equation of motion from the model. 

III. They retain only one term in the strain energy U among the additional third- and fourth- 
order terms which arise from the nonlinearities in the strain-displacement relation. 


With the above assumptions, the kinetic and potential energy functions take the form: 


T = 
U = 



(23) 

(24) 


where in accordance with the assumption II above, we have dropped the term in the strain 
energy which involves uo only. These correspond to Eqs. (39) and (45) in Ref. 7. Note, that 
the second term in (23) includes the rotary inertia contribution. Being constant, it does not 
contribute to the motion equation for w. Also, the second term in (24) is the additional term 
mentioned in III and is of third order. It represents the coupling between axial and transverse 
displacements which leads to the stiffening of the beam in bending. This term is identical to 
the third-order term in the strain energy function of Hanagud and Sarkar. 


At this point in their development, Likins et al. establish a connection between their 
approach, as we have just outlined, and “the textbook derivation for the transverse vibrations 
of beams subject to an external axial force P." As noted by them, the axial load to first 
approximation is given by: 

P = EA ^ (25) 

ox 


so that Eq. (24) can be rewritten as: 



(26) 


To proceed with the application of Hamilton’s variational principle, Likins et al. assume that 
P is time-independent and can be approximated by its steady-state value. In fact, for the 
particular problem of a beam rotating at a constant speed, the axial load P is the centrifugal 
load on the beam. Furthermore, since the latter is a known function of the prescribed ft, P 
can be calculated with: 

P(x) = \ptf(L 2 - x 2 ) (27) 


With equations (23), (26) and (27), Likins et al. derive the motion equation for the elastic 
displacement w. Due to the assumptions made in evaluating the strain energy and the axial 
load P (Eqs. (24), (25) and (27)), geometric stiffening appears as a. linear term in this equation. 

As a final comment on Ref. 7, we note that in section titled “Finite- Element Model” 
Likins et al. explicitly include a term that represents “modifications of structural stiffness due 
to spin-induced loads on the structure in its steady state (the so-called “geometric stiffness”).” 
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Let us now proceed with the developments presented by Vigneron, and Kaza and Kvater- 
nik. As we had already mentioned, both derivations are premised on expansion (21) for the 
axial displacement. As well, the kinetic energy is constructed by ignoring the effects of rotary 
inertia. This can be achieved by neglecting the “off-elastic-axis” contribution, (-yf| - z§~), 
in the expression for u and dropping the last (/*-) term in Eq. (17). Then, kinetic energy can 


be written as a function of 

u a , v, w and uj with: 


T= ^ 

J (ii] + i> 2 + w 2 — 2u 3 u j + Uy) dx 


+¥/> 

J ^r 2 + u 2 + v 2 + u 2 + 2xu , - [2 xuj] - 2u,u/) dx 



J (^—u 3 v + viij + xv + vu s — ujiij dx 

(28) 


Recall that according to (22), the foreshortening uj is a quadratic function of the transverse 
displacement gradients and, therefore, the above expression for T includes terms of third and 
fourth order. These terms are underlined in (28) with single and double lines, respectively. 
We have also singled out by enclosing in square brackets the term which is linear in uy. This 
term gives rise to the stiffening effect in the motion equation. We also emphasize that Eq. (28) 
originates from the same expression for the kinetic energy (Eq. (17)) as used by Likins et ai 
The higher order terms in it are a consequence of dividing the axial displacement u 0 into two 
components, u s and u/, and explicitly assuming a second-order form for the latter. 

Vigneron approximates T by keeping only second-order terms in (28) (including the 
term in the square brackets) as well as, setting the axial displacement u a and its derivatives 
to zero. The latter approximation corresponds to inextensibility assumption, which implies 
that the beam is modelled as axially rigid. 9 Note, that dropping the u^-terms in Vigneron’s 
formulation is not equivalent to dropping the uo~ terms in Likins et a/.’s formulation, although 
in both cases this achieves elimination of the independent axial equation of motion. Vigneron’s 
approximation yields: 

T = \p J (i) 2 + tb 2 ) dx + j (^ 2 + v 2 - [2 £i£/]) dx + pSl J (xv) dx (29) 

which, without the last term, corresponds to equation (9) of his paper. Unlike Vigneron, 
Kaza and Kvaternik do not make the inextensibility assumption and therefore keep the axial 
equation of motion in their final model. They also retain third-order terms in (28). The result 
is given by Eq. (7) in Ref. 9 which is identical to our Eq. (28) without the double-underlined 
terms. 

To determine the strain energy according to (18), Vigneron and Kaza and Kvaternik 
first evaluate the strain e xx of (19) by substituting for the axial displacement u from (21) in 
conjunction with (22). Upon dropping the third- and fourth-order terms, this produces the 
following second-order expression for the axial strain: 



where we have intentionally separated the contribution of uj. However, according to (22), we 
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have 


( 31 ) 


o*l _ i IV^y . (?m.Y 

dx 2 [V^/ \dx ) 

so that (30) is simplified to: 



or 


du 3 d 2 v d 2 w j 
= ~dx~ V d x* ~ z !h? + 2 


du 3 

dx 


d 2 v 3 2 w 
y dx 2 2 dx 2 


(33) 


The above is equivalent to Eq. (5) in Ref. 8 and Eq. (4) in Ref. 9, although neither publication 
includes a description of the intermediate step (30). It is interesting to note that £ xx evaluated 
with (32) without the second-order term, can be derived from the axial strain of (19) without 

the term, but with the “transverse” nonlinear terms. As is demonstrated by Eq. (30), 

the latter is cancelled when foreshortening uj is explicitly defined in the axial displacement u, 
which results in a linear expression for the strain. Vigneron’s strain energy can be constructed 
by using this linearized form of (32), since he retains only quadratic terms in the strain energy. 
In fact, the strain energy function given in Ref. 8 has a standard form used in the linear theory 
of elasticity. Kaza and Kvaternik formulate their strain energy with the axial strain as given 
by (32), but drop the resultant fourth-order terms. 

Because of the differences in the approximations made, Vigneron and Kaza and Kvaternik 
derive different sets of motion equations. Vigneron obtains two linear equations of motion for 
the elastic variables v and tv (Eqs. (10) and (11) in Ref. 8), the latter being identical to the 
equation derived by Likins ef al. Kaza and Kvaternik present a set of nonlinear dynamics 
equations for the three elastic deformations, u 3 ,v and w (Eqs. 8(a,b,c) in Ref. 9). They are 
nonlinear, in particular second-order, because Kaza and Kvaternik retain third-order terms in 
their kinetic and strain energies. It is important to emphasize, that the term responsible for 
the stiffening of the beam derived in Refs. 8 and 9 is a first-order term and appears via kinetic 
energy. This occurs because of the particular form assumed for the axial displacement (Eq. 
(21)) as well as, the expression adopted for the foreshortening (Eq. (22).) If, as is done by 
Likins et al., 7 one does not explicitly identify uj in the axial displacement, then the stiffening 
appears in the motion equations through the strain energy and is fundamentally a nonlinear 
term. 


To conclude this section, we draw attention to some of the observations made by Kaza and 
Kvaternik. They identify four different approaches for deriving linear or nonlinear equations of 
motion. They are: (1) “the effective applied load artifice;” (2) the use of Newton’s second law 
applied to the deformed configuration; (3) an approach in which nonlinear strain-displacement 
relations and a first-degree displacement field are used; (4) Vigneron’s approach which uses 
nonlinear strain-displacement relations and a second-degree displacement field. Kaza and 
Kvaternik show that all four appraches “make use of geometric nonlinear theory of elasticity 
either implicitly or explicitly.” They state that “for developing the equations of motion for a 
rotating beam, “the geometric nonlinear theory is necessary to obtain even the correct linear 
equations.” In their paper, Kaza and Kvaternik also discuss whether foreshortening must 
be explicitly included in the axial displacement field. They conclude that although it is not 
necessary, the approaches where the foreshortening effect is accounted for otherwise, require 
special considerations. 
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3.2 Lips and Modi, Hughes and Fung 

Lips and Modi investigate the dynamics of satellite systems composed of a central rigid body 
with flexible appendages. The base body is allowed to undergo general rotational motion and 
therefore, the dynamics model includes the rotational rigid-body equations. To illustrate the 
procedure for modelling deformation of appendages, Lips and Modi present an explicit form 
of the elastic equations for a rotating beam. In their treatment of this system, they take into 
account: (i) the offset between the attachment point of the appendage and the center of mass 
of the rigid body; (ii) variable flexural rigidity which subsumes variable modulus of elasticity 
and cross-sectional area; (iii) variable density. The basic assumptions made in deriving the 
equations for elastic displacements are similar to those made by Kaza and Kvaternik. 9 In 
particular, Lips and Modi explicitly separate the foreshortening component from what they 
refer to as the oscillation component in the assumed form of the axial displacement. As well, 
the kinetic and potential energy functions employed to construct the Lagrangian contain terms 
up to third order. 

The main subject of the work by Hughes and Fung 11 is the stability of spinning satellites 
with long flexible appendages. Therefore, they formulate the dynamics equations with a view 
to addressing this issue. The system is modelled as a spinning rigid body with appended 
beams. The rigid-body equations are formulated for small perturbations from the nominal 
spin configuration. The development of the elastic motion equations presented in Ref. 11 
is different from the previously considered works in two respects. First, Hughes and Fung 
employ deformed coordinates to describe the kinematics of the deformed beam, a fact which 
they do not state explicitly. Thus, the position of a point on the elastic axis is defined by 
[A\ Vj w] T where X = X (t ), v — u(X, t) and w = w(X y t). (Note, the corresponding undeformed 
description is [x + u,v,w] T where u = u(x,i), etc.) Another distinct feature of Hughes and 
Fung’s formulation is that they evaluate the kinetic and potential energies by integrating the 
respective appropriate integrands over the volume, which in the case of a slender beam reduces 
to the arc length s, of the deformed configuration. Thus, the energy functions are defined by 
means of the line integral /(-)ds. In general, the integration necessary to determine kinetic 
and potential energies of an elastic body can be performed over the deformed or undeformed 
configurations. However, the former is the standard choice when the kinematic description is 
given in terms of the deformed coordinates. (This representation corresponds to the Eulerian 
or spatial description of the problem.) 

Hughes and Fung incorporate the geometric nonlinearity into their formulation by ex- 
pressing the differential arc length ds with 

ds = Vi + t/ 2 + w' 2 dX (34) 

The prime in the above denotes differentiation with respect to the deformed coordinate X . We 
note that Eq. (34) is equivalent to Eq. (4) of Hanagud and Sarkar. 6 It defines the distance along 
the beam as a function of the transverse displacement gradients, when these are expressed in 
terms of the deformed coordinates. Hence, Eq. (34) is also equivalent to Eq. (1) of this paper, 
with s = x + s. To simplify the derivation, Hughes and Fung use a second-order approximation 
of ds so that ^ ^ 

/(•)dJ = jf ()dX + \J o (■) (y n + w*) dX (35) 

where L* denotes the projection of the tip on the axial coordinate axis. We also note that by 
using L in the upper limit of the integral f(-)ds, Hughes and Fung implicitly assume that the 
beam is inextensible. 
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The kinetic and potential energies derived in Ref. 11 contain second-order terms only. 
That, combined with the inextensibility assumption makes the formulation of Hughes and Fung 
similar to Vigneron’s. In fact, the kinetic energy of a single beam rotating at a nominal speed, 
deduced from Eq. (9) of their paper, is equivalent to Vigneron’s kinetic energy given by our 
Eq. (29). The form of the strain energy function is also the same. 

This completes our overview of the publications disputed in London’s and Kane et a/.’s 
commentaries. 5,l!j It is clear that these works incorporate geometric stiffening into the dynamics 
equations of a rotating beam, although through different approaches. A discussion of these 
will be given in §6 of the paper. At this point, we only note that they all account for the 
coupling between the transverse and axial deformations. Indeed, it is exactly this phenomenon 
that causes stiffening of an elastic body under certain conditions. The differences between 
the approaches lie in what we view as the mechanism for introducing the coupling effect into 
the formulation and accordingly, the stage in the derivation at which it is introduced. In 
the following two sections of the paper we discuss some of the other approaches that have 
been employed to account for the geometric stiffening in the dynamics equations of multibody 
systems. 


4 Laskin et al., Meirovitch, Banerjee et al. 


4.1 Laskin et al. 


Similar to Vigneron, Laskin et a/. 12 assume a form for the axial displacement u in which the 
displacement of points along the elastic axis, uq , is divided into two parts. We rewrite their 
Eq. (6) as: 

dv dw 


u = u, 


ijj 


+ “'-^-‘31 


(36) 


where they refer to u qa (u 0 in Ref. 12) as a quasi-steady component and u t (their v*) as a 
transient component that accounts for longitudinal vibrations. They justify this arrangement 
by arguing that it allows one to consider u< as a small, more precisely, infinitesimally small 
displacement which is of the order of the transverse displacements v and w. The quasi-steady 
component u qs may be comparatively large. 

To discretize the elastic deformation field, Laskin et al. use modal expansion similar in 
form to (3) but written for the transient axial displacement u t rather then u 0 (ui in (3)). By 
doing so, they implicitly choose u t as a generalized coordinate in their formulation. Laskin et 
al. explain this choice by saying that because u t is small, its modal coordinates can also be 
regarded as small. By constrast, if one discretizes tio = u qa + t ij, its modal expansion cannot 
adequately account for both large ( u qs ) and small (uj) components without having to include 
a large number of terms in the expansion. 

Laskin et al. employ Kane’s method to derive the motion equation of the beam. Thus, 
the next step in the formulation is to develop the generalized inertia and active forces. Starting 
with the latter, they consider two contributions — the elastic forces and the controller forces. 
In keeping with the subject of this paper, discussion is limited to Laskin et a/.’s derivation of 
the generalized elastic forces. 

Let us recall that Kane et al., and Hanagud and Sarkar derive the elastic force from the 
strain energy function. Laskin et al. take a different approach and express the elastic force 
F as a function of the stresses in the beam. This is actually a traditional way of writing the 
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elastic force contribution to the equilibrium equations of the theory of elasticity. As pointed 
out, 12 one can express F in linear elasticity with: 

F = V • cr (37) 


where cr is the stress tensor. We note that this form is also valid in nonlinear theory, provided 
the divergence operator is defined with respect to the deformed coordinates and the strains 
are small. If V is defined with respect to the undeformed coordinates, then Eq. (37) must be 
modified to take into accout the nonlinear nature of the deformation gradient. The resultant 
expresssion for F is given by Eq. (21) in Ref. 12 which we do not repeat here. In agreement 
with Kaza and Kvaternik’s comments on the Newton’s second law approach, we point out that 
this form allows for arbitrary rotations and implicitly assumes nonlinear strain-displacement 
relations. 

Since the generalized elastic forces given by (22)-(25) 12 are expressed in terms of the 
generalized coordinates, Laskin et al. must take an intermediate step of substituting for the 
stresses in their (21) from stress-strain and then strain-displacement relations, in order to 
reformulate F in terms of u q3 , u t , v and w. What is most notable about their approach is 
that, by using the “nonlinear” relation for F in terms of the stresses, they are able to derive the 
geometric stiffening component of elastic forces with a linear form of the strain-displacement 
relations. Laskin et al. emphasize that if one were to use strain energy to derive the elastic 
forces F, one would have to retain nonlinear terms in the strain-displacement relations. To add 
to their interpretation of this “paradoxical situation”, we offer this observation. The nonlinear 
terms in the strain-displacement relations give rise to third- and fourth-order terms in the 
strain energy. As was shown in §3.1, it is the third-order terms that are responsible for what is 
usually referred to as geometric stiffening. (In fact, the fourth-order terms also contribute to 
stiffening of the beam in bending.) Since the “nonlinear” formulation of F in terms of stresses 
is one-order higher then the “linear” one, it essentially provides a mechanism for incorporating 
only third-order terms in the strain energy. 

We note that geometric stiffening is represented by G tJ - and a l0 - terms in Eqs. (23)- 
(25) of Ref. 12. Although these are linear in the discrete generalized coordinates, they are 
also dependent on the spatial derivative of the quasi-steady axial displacement u q3 (see the 
definitions of G,j and a,o on p. 515 of Ref. 12, with i>o = u q*)- 

Derivation of the generalized inertia forces 12 is similar to the procedure in Ref. 1 with 
the main difference being the addition of six rigid generalized inertia forces. These appear 
because the rigid-body motion of the beam is not prescribed, but is unknown. Accordingly, 
the final equation of motion include six equations for the position and orientation of the floating 
reference frame and the differential equations for the discrete elastic coordinates. They are 
explicitly dependent on u q3 , and therefore require this axial displacement as an input. 

Laskin et al. apply their dynamics equations to a number of special cases. They specify 
u q3 as the stretch that would occur if the beam were executing its nominal or intended motion. 
This stretch is defined as a solution of a linear second-order differential equation which we 
deduce from examples discussed in Ref. 12 to be: 


EA 


d 2 u q , _ 
dt 2 


dP 

dx 


P 


(38) 


Here, P is nominal axial load on the beam and we have introduced the bymbol p to denote 
the axial load density. To generalize this scheme for the case of general rotational motions, 
Laskin et al. propose to approximate p by the axial component of the centripetal acceleration 
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(multiplied by an appropriate inertia) of a point on the elastic axis of the beam. This yields, 
in accordance with Eq. (78) in Ref. 12: 

P = p(x + U,,) (w| + (39) 

It is suggested that the above should provide a good approximation for the steady-state axial 
load in the case of rotational motions at low angular acceleration rates. Moreover, it allows 
for a closed-form solution of (38) for u qs which can then be used as an input to their dynamics 
model. 

To complete this section, we observe that Eq. (36) is analogous to (21) employed by 
Vigneron and Kaza and Kvaternik in their formulations, with the correspondence uj = — u q3 
and u s = u t . The term “quasi-steady” used by Laskin et al to qualify foreshortening can be 
interpreted to reflect the fact that this axial displacement is present even when there is no 
axial vibration, and furthermore, it exists even under static loading. Unlike what is done in 
Refs. 8 and 9, Laskin et al do not substitute for u qs in terms of v and w as in (22), nor any 
other expression. As a consequence, u q3 appears in both generalized inertia and elastic forces. 


4.2 Meirovitch 

In his 1967 book, Meirovitch 13 includes a section on the effect of axial forces in the bending 
vibration of a bar, which as he states, cannot be ignored in some cases. In this section, 
Meirovitch derives an equation of motion for the transverse displacement of the beam by 
means of extended Hamilton’s principle. Thus, expressions for kinetic energy T and work 
function W are developed. The former takes the simple form used in planar bending vibration 
problems without the axial force. Rewritten in our notation, T specified in Ref. 13 is: 

(4o) 

In evaluating the work function, he proposes to include the effect of the bending moment, the 
transverse (external) load and axial force. The first two are formulated in the same way as for 
the case without the axial force. To determine the “axial” work, the change in the horizontal 
projection of an element ds is calculated. This differential of the foreshortening is expressed 
with: 

ds — dx = 1 + (|^) dx -dx*\ (|^) dx (41) 

where the approximation results from retaining two terms in the binomial expansion of \j\ + 
Then, the work done by the axial force is: 

Wp = -*Jo P{tJ) (^) dX (42) 

It is worth to point out that in adopting the above formulation, Meirovitch makes a tacit 
assumption that the axial force is given as a known function of x and t. 

In the more recent publication, 14 Meirovitch derives a set of motion equations for a 
general flexible body in general motion. These equations are written in terms of the rigid- 
body quasi-coordinates and the continuous elastic coordinates u, v and w. Their application 
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to a system made up of a rigid hub and a beam-like flexible appendage is illustrated. He begins 
by assuming that the axial displacement can be ignored and therefore sets u = 0 apriori. The 
kinetic energy is derived in the standard manner and contains terms that are of second degree 
in the elastic variables. The strain energy includes the standard second-order contributions due 
to bending in two directions as well as the contribution due to “shortening of the projection.” 
The latter is expressed as 

"'-if [jf* f H[(®) ,+ (£) 1* (43) 

where p(x,t) is the axial component of the internal force density. We note that (43) can be 
directly compared to (42). 

Meirovitch proposes to determine the force density p from the motion equation for the 
axial displacement u, which as we had mentioned is excluded from the dynamics model. Thus, 
he defines p as a sum of the terms in the u-equation, omitting terms that involve elastic 
displacements, as well as the control force density. The resultant expression for p is given by 
Eq. (29) in Ref. 14, which we rewrite here as: 

p = p f — Vj - w 2 v 3 + o> 3 v 2 + i (u\ + w 3 )] (44) 

At this point, let us compare the above expression for the internal axial force density which is 
employed by Meirovitch to evaluate the strain energy due to the shortening of the beam with 
the corresponding expression of Laskin et al., which they use to determine the quasi-steady 
component of the axial displacement. We recall that the expression for p proposed in Ref. 12, 
as given by Eq. (39) of this paper, originates from inertial acceleration of a mass element on 
the beam’s elastic axis. In fact, it is defined as the centripetal component of the acceleration 
evaluated with u = u q> . It can be shown that Meirovitch ’s expression for p can also be derived 
from the inertial acceleration distribution, but neglecting the elastic contributions. Indeed, 
this is not surprising since a motion equation in the absence of external forces is essentially a 
linear homogeneous relation for the inertial acceleration. Thus, the “centripetal component” 
of p, given by px (u^ + in Eq. (44), is identical to p of Eq. (39) if one drops u q3 . Finally, 
by comparing (44) to (39), we observe that Laskin et al.' s approximation for the axial force 
density should apply when the translational velocity and acceleration are small relative to the 
angular velocity, in addition to the conditions stated in Ref. 12. 

4.3 Banerjee et al. 

Banerjee and Dickens 4 present a formulation for the dynamics of an arbitrary flexible body 
in large rigid-body motion. The formulation is based on the standard description of the 
deformation field where the elastic translations u, v and w are discretized by means of a modal 
expansion. The equations of motion are derived by employing the same general methodology 
as in Refs. 1 and 6. Therefore, as was done in reviewing other works, the following analysis 
focuses on the procedure to derive the stiffening terms only. 

Banerjee and Dickens introduce the notion of “motion stiffness” as a special case of the 
geometric stiffness caused by the inertia loading on the body due to its large rigid motion. As 
is noted by Banerjee and Lemak, 15 this motion-induced stiffness has its origin in the strain 
energy term, 

U NL = J e NL T a dV (45) 
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where a and e NL are 6x1 columns. One can view Uhl as the “nonlinear” part of U , since it 
arises from the nonlinear terms in the strain-displacement relations. These are represented in 
(45) by e NL , while a represents the state stress in the body. 

The key to the method proposed Ref. 4 is the observation that the stress state of the 
deformable body undergoing large translations and rotations results from the inertia loading 
on it. With that, they determine the distributions of the inertia force and torque in the body 
from the inertial linear and angular accelerations, respectively, in accordance with Newton’s 
Second Law and Euler equations. In applying this procedure, they neglect elastic terms in 
the velocity and acceleration distributions. Thus, the resultant inertia loads are expressed as 
functions of rigid velocities and accelerations. 

In the next step of their development, Banerjee and Dickens rewrite the inertia load as a 
product of a particular matrix and a column vector. For the inertia force, the matrix is 3x12 and 
is dependent on the spatial coordinates, x, y and z. The column vector contains 12 parameters 
A{, i = 1, . . ., 12, which are dependent strictly on the velocities and accelerations of the body. 
The inertia torque is factorized into a constant 3x9 matrix and a time-dependent 9x1 column 
vector of A,, i = 13,..., 21, which are also functions of rigid velocities and accelerations. 
Although not stated in Ref. 4, the motivation for this “factorization” is to separate the time- 
dependent component of the inertia forces from the space-dependent or constant part. By 
doing so, Banerjee and Dickens are able to construct the geometric stiffness term in the motion 
equations in two stages. The first stage produces 21 geometric stiffness matrices denoted by 
5(*M These can be assembled with the standard finite-element procedure from the constant or 
space-dependent matrix factors of the inertia loads, prior to evaluation of the motion equations. 
The second stage involves calculating the geometric stiffness term from A, and the discrete 
elastic coordinates. 

Let us now comment on the relationship between the approach of Banerjee and Dickens 
and the other methods. First, we observe that the “nonlinear” strain energy of Eq. (45) is 
“exact” in the context of small strain deformation. Moreover, it can be rewritten in terms 
of displacement variables if one expresses stresses in terms of strains and then substitutes for 
strains in terms of displacements. Following this procedure, one will obtain Uhl in the form 
of third- and fourth-order terms in the displacement gradients. These were mentioned in our 
discussion of Hanagud and Sarkar’s work. 

It can be shown that Banerjee and Dickens’ expression for the axial inertia force is 
equivalent to Meirovitch’s axial component of the internal force density. The latter is given by 
our Eq. (44) and the former can be obtained from Eqs. (28) and (29) of Ref. 4. Substituting 
for xj — x, X 2 = y = 0, 13 = z = 0 and u; = v,, 1^+3 = w, for i = 1, . . . , 3, the axial inertia 
force /,* of Banerjee and Dickens takes the form: 

ft = -[Ai+xA 4 ]dm 

= - [vj + u-iv-j, - (j 3 v 2 - 1 (u>l + u>|)| dm (46) 

Since Meirovitch’s internal force density must be interpreted as the internal force per unit 
length, the equivalence between (46) and (44) follows if one substitutes for dm with p dx. 
Then, we have 

ft = pdx (47) 

where p is as defined by Meirovitch (Eq. (44)). Furthermore, the axial stress in the beam is 

<7 rr (x) = \ j P dx ( 48 ) 
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( 49 ) 


and if the nonlinear part of the axial strain (19) is simplified to 


fjVL.n 



Eq. (45) applied to the beam reduces to: 


U NL 


£nL,xx&x 


J 


(50) 


The above, which we derived from Unl of Ref- 4 > s identical to Meirovitch’s strain energy due 
to shortening of the projection. In this light we propose to interpret Meirovitch’s formulation 
of the strain energy due to shortening of an elastic beam as a particular application of Banerjee 
and Dickens’ formulation. 

To conclude this section, we would like to draw attention to a fundamental approximation 
made in the development of Banerjee and Dickens. It is that the stresses which contribute to 
stiffening of the body are only those that are caused by the inertia loading. We believe that the 
stress in (45) must represent the complete stress state in the body which arises from the total 
loads— inertia and external. Indeed, only then can Eq. (45) be consistent with the general 
formulation of the strain energy. The approximation for a suggested in Ref. 4 is certainly 
valid in the absence of external forces on the body. However, it may not be appropriate for 
multibody systems in which each body is acted upon by the “externaT (to it) interbody forces, 
even in the absence of external forces on the whole system. 


5 Discussion 

Based on our review, we propose classification of the different approaches to model geometric 
stiffening according to these two criteria: (i) kinematic description of the deformation field 
(criterion DF); (ii) the formulation of the strain energy function (criterion SE). We feel that 
these represent two most general criteria that can fully characterize a particular approach. The 
first one is related to the assumed displacement field, which in turn determines the generalized 
coordinates employed in deriving the motion equations. The second criterion defines the form 
of the elastic force term in the motion equations. 

According to each criterion, we have identified the following possible cases. 

Criterion DF. It is suggested that the deformation field can be described by any of the 
following sets of variables: 

DF-1: Three independent elastic translations {u,u,w}. This description is employed by 
Likins et al., Hanagud and Sarkar, Meirovitch, and Banerjee et al. 

DF-2: Three independent elastic translations {u^w, tn}, where recall u , is the axial dis- 
placement which results strictly from the stretching of the deformable body. This 
description is used by Vigneron, Kaza and Kvaternik, Lips and Modi, and Laskin 
et al. 
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DF-3: Three independent but non-orthogonal elastic deformations {s,u,tn} where s de- 
notes the stretch. These coordinates are employed by Kane et al. 

Criterion SE. We distinguish two basic formulations of the strain energy: 

SE-1: The total strain energy is formulated as a function of strains (or displacements) 
only, where strains are defined with the nonline ar strain-displacement relations. 
This formulation is employed by Likins et al., Vigneron, Kaza and Kvaternik, Lips 
and Modi and Hanagud and Sarkar as well as, Hughes and Fung, and Kane et al. 
It can be viewed as a displacement formulation. 

SE-2: The strain energy is subdivided into a “standard linear” contribution and a “nonlin- 
ear” part, where the latter is constructed from the stresses or forces in the body and 
the nonlinear part of strain. This formulation is used by Meirovitch, and Banerjee 
et al. and can be considered as a force formulation. 

Note that Laskin et al. do not compute the strain energy, but derive the elastic forces 
directly from the stress state. However, since they eventually reformulate these in terms 
of strains, their method is fundamentally a displacement formulation. 

Let us comment on the three options for the kinematic description of the deformation 

field. 


DF-1. The set {u,u,u;} is the standard set of elastic displacements used in both linear and 
nonlinear theories of elasticity as well as, structural analysis. It requires no apriori 
decisions on the foreshortening part of u and in that sense is general. 

Since u, v and w represent deformations along three orthogonal axes, they clearly rep- 
resent independent degrees of freedom. However, contrary to what has been previously 
stated, that does not preclude the fact that part of u, in particular, the foreshortening 
component, is dependent on v and w. The total axial displacement will remain an in- 
dependent variable as long as it includes the axial displacement tt 4 . Furthermore, it is 
not necessary to either separate u into these two components, or to explicitly account 
for the coupling between the axial and transverse displacements. In addition, we believe 
that discretizing u with a standard expansion does not imply “premature linearization.” 
To support this point, let us consider the foreshortening part of u, which for the present 
purpose we simplify to: 



If the transverse deformation v is expanded as in (3), then upon substitution, uj becomes: 


u j = \ E £ 

j=l k = 1 


jTW)W) 


(52) 


The above clearly represents a summation of terms <f>/,i(x) qj,i(t) where qf,i(t) = qj(t) qk( 0 
and the space-dependent basis functions are defined accordingly. This summation is of 
the same form as the standard expansion for u. 

As exemplified by the formulations of Likins et al. and Hanagud and Sarkar, geometric 
stiffening can be modelled with this option, provided one incorporates the nonlinear 
strain-displacement relations in the evaluation of strain energy and, hence, the elastic 
forces. 
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DF-2. In this case, the foreshortening component of u is explicitly separated from u, in the 
axial displacement field. It may be specified in terms of v and w or left as a parameter 
to be defined by the user. Either case, however, involves making an approximation for 
u/, although in the first option it is “known” and quantified prior to deriving the motion 
equations. Recall that Laskin et al. have argued in favor of this option based on their 
claim that displacements u s , v and w are of the same order, while u/ is comparatively 
large. In response to this statement, we would like to draw attention to the results 
presented in Ref. 6, where Hanagud and Sarkar display the time- histories of axial and 
transverse displacements from a simulation of a beam spin-up problem. As can be seen 
from Figures 2a and 2b in Ref. 6, the axial deformation «o = «, — «/ is almost two orders 
of magnitude smaller than the transverse displacement v. Therefore, u, cannot be of the 
same order as v, nor can uj be significantly larger than v. 

DF-3. Employment of the stretch variable is unconventional, and certainly is not a common 
choice, if made at all. The main advantage of using this variable instead of u is that 
the strain energy retains its (“linear”) quadratic form. However, the resulting expression 
for kinetic energy (or generalized inertia forces) is more complicated than that based on 
DF-1 and DF-2 descriptions. 

As was demonstrated in the previous sections (at least we hope it was), it is a particular 
combination of the coordinates and the strain energy formulation that determines how the 
geometric stiffening is incorporated into the motion equations as well as, what form it appears 
in. Thus, we will now comment on the two strain energy formulations taken in combination 
with the different DF options and point out some of the advantages and/or disadvantages of 
the resulting approaches. 

The main advantage of the SE-1 formulation of the strain energy for any description 
of the deformation field is that one is not required to make any additional assumptions or 
approximations. 

SE-l/DF-1. In this approach, taken by Hanagud and Sarkar, the stiffening term appears in the motion 
equations through the strain energy and is a nonlinear function of elastic coordinates. In 
particular, it can be factored into an (elastic)coordinate-dependent second-order stiffness 
matrix and a column vector of elastic coordinates. 

As pointed out by Banerjee and Lemak, evaluation of this stiffness term requires that 
the stiffness matrix be updated at each time step in the simulation, which may be com- 
putationally costly. However, the update does not involve iterations, but is a simple 
functional evaluation. It should also be pointed out that this approach requires that the 
axial elastic equation be included in the dynamics model. This is likely to have an ad- 
verse effect on the computational efficiency of the simulation, since the axial deformation 
is usually a high-frequency component. 

Finally, we note that the present method corresponds to the third approach identified 
by Kaza and Kvaternik and as they comment is “the one usually employed for a general 
three-dimensional rotating body.” However, contrary to Kaza and Kvaternik’s conclu- 
sion, we believe that it does not require special consideration, but is the most general. 
Moreover, this approach can be employed to extend the existing “small” deflection dy- 
namics formulations to incorporate “large” deflection theories. 

SE-l/DF-2. With this approach, the foreshortening term is always present in the kinetic energy 
expression and may or may not appear in the strain energy, depending on the additional 
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approximations made. In the simplest case, as in Vigneron’s formulation, the resultant 
stiffening term takes a linear form. 

This method of determining the stiffening term allows one the option of dropping the 
axial deformation from the model— not a trivial advantage from the computational point 
of view. However, it may be less accurate than the previous method, because of the 
approximation made in assuming a form for the foreshortening displacement. 

SE-l/DF-3. This approach produces a linear geometric stiffening term via kinetic energy (or gen- 
eralized inertia force). Compared to the SE-l/DF-1 option, it can be just as accurate, 
but not as general since the stretch variable can be defined only for a particular type of 
elastic bodies. 

As was shown in §4.3, the SE-2 formulation of the strain energy requires an approximation 
for the stress field in the body, and hence, in contrast to SE-1, is inherently approximate. 
Among the works presented in this paper, this formulation has only been used with the DF-1 
description of the deformation field. In this case, the geometric stiffening term results from the 
“nonlinear” strain energy, as in Refs. 4, 15 and 14. With appropriate assumptions, it is linear 
in elastic variables, but involves rigid-body accelerations and velocities. As for its efficiency, 
the algorithm presented in Ref. 14 also requires an update of the geometric stiffness matrix (see 
Eq. (37) in Ref. 4), because of the time-dependent quantities A{. Moreover, the final motion 
equations no longer have a symmetric coefficient matrix (because of acceleration dependency of 
the stiffening term), thereby making evaluation of accelerations computationally more costly. 

To summarize, we believe that a description of the deformation field in terms of u, v, w 
combined with the displacement formulation of the strain energy is the most accurate and 
general approach. It does not require an approximation of foreshortening, nor the stress 
state of the body — two critical advantages for applications to multibody systems. A definitive 
statement on the efficiency of this approach and how it compares to, for instance, Banerjee and 
Dickens’ procedure can only be made through implementation of both methods. Moreover, we 
would expect the relative efficiency of the two formulations to vary depending on the complexity 
of the system and the number of elastic degrees of freedom in the model. 

6 Concluding Remarks 

In this paper, we presented an exposition of several approaches to model the geometric stiff- 
ening effect for dynamics simulation of flexible-body systems. Our review included 11 papers 
published in the period from 1973 to 1991. Although it does not represent a complete litera- 
ture review of the works that have addressed this issue, it covers a wide range of formulations 
developed for the problem. 

In reviewing these works, we have idenitified two key characteristics of the different meth- 
ods which allowed us to put forward a general classification for them. We have also established 
the interrelationships between the various approaches, provided a number of clarifications and 
new interpretations and offered our opinions on their benefits. It is hoped that this work will 
contribute to a better understanding of the origin of geometric stiffening and how this effect 
can be incorporated into the dynamics model of a flexible body undergoing large rigid-body 
motion. 
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